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ABSTRACT 


A  deephole  array  consisting  of  12  short -period  triaxial  seismometers  was  operated 
at  a  test  site  near  Grapevine,  Texas  (GVTX) ,  until  23  February  1968.  The  system 
was  then  put  in  a  standby  status  until  additional  software  could  be  developed 
to  provide  additional  filtering  capabilities.  The  complete  system  was  shut 
down  and  all  of  the  instrumentation  and  equipment  was  moved  into  storage  areas 
at  the  Geotech  plant  in  Garland,  Texas,  during  the  last  week  of  May. 

The  information  gathered  from  the  short -period  triaxia]  array  at  the  GVTX  site 
was  used  to  study  short -period  noise  and  signals,  to  develop  signal  enhance¬ 
ment,  and  filtering  techniques.  Most  of  the  time  was  used  in  trying  techniques 
‘'used  on  least -mean-square  filtering.  In  general,  it  can  be  concluded  that 
because  the  noise  is  not  stationary,  these  filters  degrade  too  rapidly  to  be 
of  practical  interest  for  online  processing.  Several  non-optimum  filtering 
techniques  were  also  tried  and  were  found  to  be  as  effective  as  the  optimum 
filters. 
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FINAL  REPORT,  PI  ,JECT  VT/7703 
DEEP  WELL  ARRAY  OPERATIONS 


1.  INTRODUCTION 


This  report  discusses  the  work  accomplished  and  the  problems  encountered  during 
the  field  operations  program  conducted  under  Project  VT/7703.  The  work  reported 
herein  covers  field  operation  of  a  12-component  short-period  triaxial  array 
system  (4  vertical  components),  modifications,  improvements,  programming,  signal 
and  noise  analysis,  and  processing  from  1  March  1967  through  31  May  1968.  It 
is  submitted  in  compliance  with  paragraph  5,  Reports  of  the  Statement  of  Work 
to  be  Done,  Project  VT/7703.  This  project  was  under  the  technical  direction 
of  the  Air  Force  Technical  Applications  Center  (AFTAC)  and  under  the  overall 
direction  of  the  Advanced  Research  Projects  Agency  (ARPA) . 

This  report  is  generally  presented  in  the  same  sequence  as  the  tasks  in  the 
Statement  of  Work.  A  copy  of  the  Statement  of  Work  is  included  as  Appendix  1. 


2.  FIELD  OPERATIONS  PROGRAM 


2.1  OPERATION  OF  TRIAXIAL  ARRAY 

The  short-period  triaxial  deephole  array  and  the  associated  equipment  that  was 
assembled  for  Project  VT/5051,  was  moved  from  the  site  near  Apache,  Oklahoma 
(APOK),  to  the  site  near  Grapevine,  Texas  (GVTX) ,  see  figure  1,  during  the 
first  two  weeks  of  March,  1967.  All  of  the  downhole  instruments  were  disassembled, 
thoroughly  checked,  and  marginal  components  replaced  where  necessary  oefore 
installation.  The  downhole  instrumentation  consists  of  four  short-period  triax 
seismometers  that  have  three  sensitive  modules  directionally  orientated  120 
degrees  to  each  other.  Figure  2  is  a  photograph  of  one  of  these  modules  and 
figure  3  shows  an  outline  drawing  of  one  seismometer  assembly.  The  instruments 
were  operated  in  the  GVTX  deephole  at  the  following  depths: 

DH4  Surface 
DH3  1067  m 

DH2  1982  m 

DH1  2915  m 

With  the  exception  of  moisture  getting  into  a  downhole  connector,  all  of  the 

instruments  were  checked  out  and  calibrated  without  incident.  Satisfactory 
operation  of  all  downhole  instruments  continued  until  the  middle  of  July  when 
the  bottom  seismometer  in  the  deephole  started  to  leak.  This  required 
pulling  all  of  the  deephole  instruments.  The  leak  was  found  in  the  switching 
module.  While  it  was  being  repaired,  a  routine  inspection  revealed  a  heat 
distorted  calibration  coil  in  module  No.  1  of  the  bottom  seismometer.  After 
repair  and  replacement,  the  array  was  lowered  to  the  desired  depths  and  put 
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Figure  3.  Triaxial  seismometer  outline  drawing 


ii'to  operation.  During  calibration,  the  mass  position  motor  in  No.  1  module 
of  seismometer  No.  2  failed.  This  required  removing  seismometer  No.  3  and 
seismometer  No.  2  from  the  deephole.  The  motor,  and  Bendix  flexures  damaged 
when  removing  seismometer  No.  2  from  the  well,  were  replaced.  Check  out  after 
installation  revealed  that  the  cable  attached  to  No.  3  seismometer  showed  low 
insulation  resistance  from  the  electrical  conductors  to  ground.  This  required 
No.  3  seismometer  to  be  removed  from  the  well  for  cable  repair.  A  cable  service¬ 
man  was  called  in  to  repair  this  cable  and  another  known  short  in  a  spare  cable 
on  a  winch  truck  not  being  used  at  that  time.  Number  3  seismometer  was  then 
attached  to  the  repaired  cable  on  the  winch  truck,  previously  not  in  use,  and 
lowered  to  desired  depth.  All  instruments  were  operational  the  last  week  of 
August,  1967,  and  have  continued  to  operate  satisfactorily  until  23  February  1968 
when  the  entire  system  was  put  in  a  standby  status. 

All  of  the  peripheral  analog  equipment,  including  the  film  recording  equipment 
operated  satisfactorily.  The  failures  that  did  occur  were  caused  mainly  by 
wear  from  normal  operation.  These  failures  were  anticipated  and  spare  parts 
were  available,  thereby  keeping  downtime  on  these  pieces  of  equipment  to  an 
absolute  minimum. 

After  digitizing  several  magnetic  tapes  of  data  for  noise  studies,  a  decision 
was  made  to  increase  the  gain  of  the  system.  This  was  accomplished  by  decreasing 
the  input  impedance  of  the  multiplexer  in  the  computer.  After  digitizing 
the  required  amount  of  data,  the  gain  was  raised  again.  This  time,  operational 
amplitifers  (op-amps)  were  added  to  the  system  between  the  phototube  amplifiers 
and  the  data  control  modules  as  shown  in  figure  4.  The  multiplexer  inputs  were  then 
restored  to  their  original  configuration  with  the  system  gain  increased  about 
35  dB.  The  first  week  of  October,  the  initial  op-amps  were  removed  from  the 
system  and  were  replaced  by  two  op-amps  in  series  in  each  data  channel.  The 
total  additional  gain  was  now  60  dB.  Signal -to-noise  ratios  for  the  system 
with  the  op-amps  in  operation  are  listed  in  table  1. 


2.2  COMPUTER  OPERATION  AND  PROBLEMS 

The  majority  of  the  system  downtime  can  be  attributed  to  the  problems  that  have 
existed  in  the  Ambilog  200  computer.  During  the  4-month  period  from  the  third 
week  of  March  1967,  to  1  August  1967,  the  computer  was  considered  operational 
337  hours  as  far  as  listing  and  editing  programs  was  concerned.  The  computer 
was  capable  of  digitizing  raw  triaxial  data  and  performing  system  response 
calibrations  for  approximately  116  hours  of  the  337  operational  hours  mentioned. 

It  was  not  possible  to  digitize  data  for  the  full  116  hours  as  software  programs 
had  to  be  checked  and  debugged,  and  system  calibrations  had  to  be  made  before 
data  could  be  digitized. 

The  major  problems  resulting  in  computer  downtime  follows. 

1.  21  March  1967  Analog  to  digitial  converter  (ADC)  failed,  oscillations 

on  some  channels,  others  putting  out  random  spikes. 

The  ADC  section  was  sent  to  the  manufacturer  in 
Fullerton,  California,  for  repair  and  calibration. 

The  Fullerton,  California,  division  was  unable  to 
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Figure  4.  Triaxial  array  system,  block  diagram 


repair  the  unit.  They  sent  the  complete  section 
to  the  factory  in  Boston,  Massachusetts,  where  it 
remained  for  approximately  6  weeks  before  repairs 
were  completed.  The  unit  was  received  11  May  1967; 
after  two  days  of  testing,  problems  still  existed  in 
the  ADC  unit  and  it  was  sent  back  to  the  manufacturer 
again  for  repair.  It  was  received  from  the  factory 
on  12  June,  checked  out,  and  this  time  performed 
satisfactorily. 

2.  14  June  1967  Relay  failure  in  the  Dura  typewriter  caused  semi¬ 

conductor  failures  in  the  OP I  section  of  the  computer. 

The  replacement  components  were  received  and  installed 
and  the  computer  was  operational  again  on  20  June. 

3.  21  June  1967  A  microswitch  in  the  Dura  typewriter  failed  and  again 

caused  semiconductor  failure  in  the  OPI  section  of 
the  computer.  Repairs  were  made  and  the  OPI  section 
readjusted.  The  computer  was  operational  again  on 
23  June. 

4.  3  July  1967  Another  microswitch  failure  caused  the  computer  to 

stick  in  the  "type-out"  stage.  A  new  microswitch 
was  installed  and  the  computer  was  again  operational 
on  the  afternoon  of  3  July. 

5.  7  July  1967  Another  Dura  typewriter  failure  occured  and  this 

time  the  Dura  serviceman  removed  the  old  typewriter 
and  installed  a  rebuilt  typewriter.  After  the  type¬ 
writer  exchange  was  made  and  checked  out,  it  became 
obvious  that  another  problem  existed.  This  problem 
was  found  in  the  voltage  reference  section.  Two 
shorted  transistors  were  found  on  the  Ref  8  card. 

After  repair  and  adjustments  were  complete,  the 
computer  was  operational  on  17  July. 

6.  26  July  1967  This  failure  was  in  the  SCU  section.  The  serviceman 

was  called  and  arrived  on  site  on  28  July  and  replaced 
a  shorted  transistor.  The  computer  was  operational 
on  29  July. 

7.  15-23  August  1967  The  Adage  serviceman  performed  complete  ambilogical 

calibrations.  Several  marginal  components  were 
discovered  and  replaced  while  performing  the  calibration. 
Calibration  was  delayed  for  2  days  because  an  0P4 
printed  circuit  board  had  to  be  shipped  from  the 
factory  in  Boston,  Massachusetts,  to  replace  one  that 
would  not  meet  specifications. 

8.  20  February  1968  Intermittant  operation  was  detected  in  the  accumulator 

section.  The  trouble  was  found  to  be  a  cold  solder 
joint  on  a  printed  circuit  board.  The  computer  was 
operational  the  afternoon  of  20  February  after  4  hours 
of  downtime. 
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Table  1.  Signal-to-noise  ratios  for  deephcle 
_ triax  array  operating  at  GVTX _ 


Seismometer 

Module 

Seismic 

background 

(counts) 

Dummy 

PTA 

(counts) 

S/N, 

background 
dummy  PTA 

1 

1 

3950 

560 

7.1:1 

2 

5500 

570 

9.6:1 

3 

3800 

500 

7.6:1 

2 

4 

6300 

590 

11.1:1 

5 

4250 

300 

14.1:1 

6 

4200 

550 

7.6:1 

3 

7 

6100 

500 

12.0:1 

8 

8000 

450 

17.5:1 

9 

8750 

600 

14.5:1 

4 

10 

10100 

650 

15.5:1 

11 

7750 

350 

22.0:1 

12 

9450 

250 

38.0:1 

These  tests 

were  made 

at  operational  magnifications  which 

recorded  a 

normal  seismic  background  level  for  each  triax  seismograph,  in  counts 
(plus  to  minus)  for  a  15-bit  binary  digital  form. 


Computer  downtime  had  become  such  a  serious  problem  by  the  middle  of  June  that 
it  became  necessary  to  initiate  a  3-month  maintenance  contract  with  the  computer 
manufacturer.  This  arrangement  worked  so  well  that  the  contract  was  kept  in 
effect  for  a  total  of  nine  (9)  months.  Reliable  and  routine  operation  of  the 
computer  was  possible  after  23  August  when  the  ambilogical  sections  were  calibrated. 

Several  engineering  changes  have  been  made  to  the  computer  during  this  contract 
period.  Seven  (7)  engineering  changes  have  been  made  by  Adage,  the  computer 
manufacturer.  One  (1)  engineering  change  was  made  by  Fabritek,  the  manufacturer 
of  the  computer  memory.  Two  (2)  corrections  were  made  to  the  Potter  tape  deck 
by  the  Adage  serviceman.  Twelve  (12)  other  corrections  of  problems  were  made 
by  the  Adage  serviceman  and  Geotech  personnel  during  normal  computer  downtime. 

These  corrections  and  adjustments  did  not  result  in  any  computer  downtime 
during  normal  working  hours.  Figures  5,  6,  and  7  show  the  rack  configuration, 
the  computer  console,  and  the  tape  deck  as  installed  in  the  instrumentation  van. 


2.3  FILM  AND  DIGITAL  TAPE  RECORDING 

Considerable  signal  and  noise  data  have  been  digitized  and  recorded  on  film  and 
tape  during  the  operational  period  of  this  project.  Some  of  the  first  digitized 
data  recorded  on  tape  was  of  low  amplitude  and  had  low  level  spikes  in  it  caused 
by  the  tape  advance  mechanism  in  the  tape  deck.  This  particular  problem,  and 
a  similar  problem  of  intermittaat  spikes,  wore  corrected  and  all  of  the  data 
since  the  last  of  August,  1967,  is  free  from  detrimental  responses. 
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COMPUTER  RACKS  SEISMOGRAPH  RACKS 


2.4  OPERATE  THE  ONLINE  COMPUTER  TO  TEST  AND  COMPARE  REAL  TIME  VERTICAL  ARRAY 
SIGNAL  ENHANCEMENT  DATA  PROCESSING  TECHNIQUES 

The  hardware  problems  encountered  during  the  early  part  of  the  project  delayed 
the  processing  of  data  with  existing  programs;  however,  a  general  filter  program 
was  written  and  ehecked  out.  This  program  will  add  any  set  of  filter  weights 
to  10  channels  of  seismic  data  and  output  the  summation  of  the  filtered  data 
and  the  raw  data  on  the  film.  The  type  of  filtering  performed  depends  on  the 
method  by  which  the  filter  weights  are  computed.  Initial  tests  were  run  using 
Weiner  optimum  filter  weights. 

The  online  filter  program  was  checked  out  and  several  days  of  data  were  recorded 
on  film.  Figure  8  shows  the  operational  results  of  this  program.  Traces  5, 

6,  and  7,  counting  from  the  top,  are  identical  traces  to  1,  2,  and  3,  respectively, 
except  that  they  have  been  digitized  prior  to  online  processing.  They  are  then 
discriminated  and  recorded  for  comparison  with  the  raw  analog  data.  These  are 
the  results  of  processing  three  channels  of  data  using  149  weights  per  channel. 

Only  the  three  deepest  seismometer  outputs  were  used  because  the  high  noise 
levels  close  to  the  surface  degraded  the  output  when  included. 

In  order  to  permit  simultaneous  recording  of  data  on  film  and  magnetic  tape, 
and  to  incorporate  the  additional  filter  routine  requested  by  the  Project 
Officer,  the  above  filter  program  was  rewritten.  The  program  had  not  been 
^completely  checked  out  at  the  contract  termination  date.  Data  was  not  being 
properly  recorded  on  magnetic  tape  even  though  the  magnetic  tape  read-write 
routines  were  operating  properly. 


3.  ANALYSIS  OF  TRI AXIAL  ARRAY  DATA 


Analysis  of  triaxial  data  does  not  involve  any  new  principles  compared  to  tht. 
vertical  component  deephole  arrays.  The  signal  waveforms  change  with  depth  on 
the  horizontals  as  well  as  on  the  verticals.  The  change  of  signal  waveforms 
from  level  to  level  and  component  to  component  can  again  be  described  in  terms 
of  relative  simple  filters. 

The  noise  on  the  horizontals  is  different  from  the  noise  on  the  verticals;  and 
contributions  from  both  Rayleigh  and  Love  modes  can  be  expected  on  the  hori¬ 
zontals.  As  a  result  of  angular  spread  in  the  directions  of  approach  of 
surface  waves,  the  coherence  between  che  vertical  and  horizontal  components 
and  the  two  perpendicular  horizontals  is  expected  to  be  very  low. 

The  directional  orientation  of  the  triaxial  seismometers  was  determined  using 
several  large  teleseisms.  This  determination  seems  to  be  fairly  accurate,  the 
deviations  of  direction  is  less  than  10  degrees. 


3.1  DETERMINATION  OF  THE  GHOSTING  FILTERS  FOR  HORIZONTAL  COMPONENTS 

Determination  of  the  effect  of  ghosting  on  the  signal  (P  wave)  waveforms  at 
depth  can  be  done  similarly  to  the  methods  used  for  vertical  component  arrays. 
There  are,  however,  some  complicating  factors.  While  the  waveforms  on  the 
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verticals  can  be  approximated  well  by  considering  only  incident  and  surface 
reflected  P  waves,  for  the  horizontals,  the  surface  reflected  SV  wave  contributes 
significantly  to  the  shape  of  waveforms.  The  arrival  time  of  the  SV  wave  and 
its  amplitude  depends  critically  on  the  Poisson  ratio  (or  S  wave  velocity)  of 
the  medium  and  the  angle  of  incidence  of  the  wave.  No  direct  measurement  of  S 
velocity  is  available  at  the  present  deephole  sites.  In  order  to  determine  the 
average  Poisson  ratio  for  a  site,  the  following  procedure  was  used: 

A  set  of  filters  was  derived  to  simulate  the  effect  of  ghosting  at  the  deepest 
radial  component  by  filtering  the  near  surface  vertical  component  recording  of 
a  seismic  event.  Various  values  of  Poisson  ratios,  o,  were  used  to  derive 
the  individual  filters  in  the  set.  The  resulting  set  of  simulated  deephole 
radial  components  was  compared  to  the  actual  deephole  radial.  The  filter  designed 
with  the  "best"  average  Poisson  ratio  gives  the  waveform  which  agrees  most 
favorably  with  the  actual  deephole  radial  components.  Such  an  analysis  was 
performed  for  both  the  Grapevine  and  Apache  sites.  Figure  9  shows  the  output 
of  several  filters  with  varying  Poisson  ratios  obtained  for  the  Greeley  nuclear 
explosion  as  recorded  at  Apache,  Oklahoma.  It  shows  that  the  best  agreement 
found  is  o«  0.3.  Fairly  good  agreement  is  found  also  for  o=  0.35.  Vertical 
waveforms  (not  shown)  found  for  this  Poisson  ratio,  however,  do  not  agree  with 
the  observed  waveforms.  For  Grapevine,  the  procedure  yielded  the  same  value 
O  =  0.3. 

To  test  the  adequacy  of  filters  to  simulate  ghosting,  a  set  of  deephole  traces 
were  derived  from  the  near  surface  vertical  using  the  Poisson  ratio  o  =  0.3. 

Figure  10  shows  the  rrsults  for  the  Greeley  explosion  recorded  at  Apache, 

Oklahoma.  The  figure  shows  that  the  waveforms  agree  very  well.  For  this 
case,  an  elastic  halfspace  model  was  used  with  only  incident  P,  reflected  P, 
and  SV  considered. 

Figure  11  shows  the  results  of  a  similar  test  for  the  Grapevine,  Texas,  site 
for  a  Nicaragua  earthquake.  Again,  good  agreement  in  the  waveforms  is  found. 

The  filters  used  for  this  test  were  cauculated  by  computing  the  response  of 
a  layered  elastic  medium  by  Haskell's  method.  A  computer  program  written  by 
Hannon  was  modified  to  yield  the  response  as  a  function  of  depth  (Hannon,  1967) . 
The  frequency  response  functions  were  inverse  Fourier  transformed  and  convolved 
with  the  surface  vertical  trace. 

In  sr  eral  previous  processing  schemes,  deghosting  was  applied  to  deephole 
trav  s  (see  final  Report,  Project  VT/5051,  TR67-3,  1967).  Deghosting  filters 
for  vertical  components  turned  out  to  be  especially  simple  two-point  recursive 
filters.  This  is  not  true  for  the  radial  components.  The  ghosting  filters 
for  radial  components  (filters  which  derive  the  ghosted  radial  trace  from  the 
incident  P  waveform)  are  not  minimum  phase  and  inverse  filtering  must  be 
implemented  by  forward  and  reverse  recursion. 

Although  all  waveforms  can  be  reduced  to  an  impulse  through  filtering  recursively 
by  the  waveform  itself,  in  general,  recursive  filtering  leads  to  instability. 

This  means  that  the  impulse  response  of  the  recursive  filter  increases  without 
bounds.  Any  small  deviation  from  the  exact  waveform  can  initiate  this  unstable 
condition.  Recursive  filters  are  stable  only  if  the  filter  used  is  minimum 
phase  or  minimum  delay,  and  are  unstable  if  they  are  maximum  or  mixed  delay. 
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Figure  9.  Determination  of  Poisson  ratio 
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Original  and  derived  traces  for  the  Greeley  event,  APOK 


TRACES  DERIVED 
FROM  SURFACE  VERTICAL 


Figure  11.  Original  and  derived  traces  (Haskell's  method)  Nicaragua 
event,  GVTX 
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It  is  easy  to  determine  to  which  group  a  particular  filter  belongs.  If  the 
filter  weights  are  slq,  a^,  a2  ■  .  .  and  the  polynomial  (Z  transform)  is 

P(Z)  -  aQ  +  aj  Z  +  a2Z^  +  .  .  .  +  a^Zn 

and  the  roots  of  this  polynomial  are  r^  .  .  .  rn  (real  or  complex  conjugate 
pair).  The  filter  is  minimum  delay  if  [r^j  >  1  for  all  roots.  The  filter  is 
maximum  delay  if  [r^]  <  1  for  all  roots. 

Mixed  delay  has  roots  both  larger  and  smaller  than  zero  in  absolute  value. 

Our  filters  are  mixed  delay,  therefore,  if  applied  recursively,  they  are 
unstable  and  the  filter  polynomial  can  be  factored  into  two  polynomials, 

P (Z)  =  PX(Z)  P2(Z) 

where  P^(Z)  is  minimum  delay  and  P2(Z)  is  maximum  delay. 

P,(Z)  can  be  applied  recursively  without  any  difficulty;  P2(Z)  cannot;  but 
it  can  be  used  in  reverse  recursion  if  we  replace  Z  by  1/Z  which  is  equivalent 
to  reversing  the  time  scale  and  thus  P2 (1/Z)  becomes  minimum  phase. 

A  detailed  discussion  of  these  concepts  may  be  found  in  the  references.  For 
our  case,  the  forward  and  reverse  recursive  filters  are  derived  by  taking 
the  root  of  the  Z  transform  polynomial  where  ag,  a^  .  .  .  an  are  the  filter 
weights  for  the  filter  which  derives  the  particular  deephole  trace  from  the 
surface  vertical.  The  roots  are  separated  into  two  groups,  one  with  absolute 
values  smaller  tnan  unity,  | r 1 |  <  1,  and  the  other  with  absolute  values  larger 
than  unity  |rj|  >  1,  forming  two  (Z)  polynomials; 

PX(Z)  =  (Z  -  rx)  (Z  -  r2)  .  .  .  (Z  -  rR) 

and,  P2(Z)  =  (Z  -  rj)  (Z  -  r^)  .  .  .  (Z  -  rm) ,  m  +  k  =  n. 

The  coefficients  of  these  polynomials  are  the  desired  filter  weights.  After 
normalization,  the  coefficients  of  the  first  polynomial  are  used  as  inverse 
recursive  filters  taken  in  inverse  order  and  those  of  the  second  polynomial 
are  used  as  forward  recursive  filters.  A  computer  program  for  factoring  Z 
transforms  in  this  manner  has  been  written.  It  was  found  that  many  roots  of 
the  Z  polynomials  for  horizontal  ghosting  filters  are  extremely  close  to  the 
unit  circle  in  the  complex  plane.  This  indicates  that  the  inverse  recursive 
filters  would  be  very  close  to  instability. 

Figure  12  shows  some  deghosted  traces  using  the  Greeley  recordings  at  APOK. 

No  difficulty  was  encountered  in  deghosting  the  deephole  vertical  traces.  This 
is  not  the  case  for  the  horizontals.  For  the  deepest  radial  component,  the 
forward  and  reverse  recursion  filters  were  found  to  ring  severely,  as  indicated 
by  the  location  of  roots,  and  the  best  results  were  obtained  by  using  a  150 
point  Wiener  filter.  The  waveforms  obtained  from  the  horizontals  are  poor. 

It  is  desirable  to  deghost  the  deephole  traces  at  some  point  during  the 
processing  sequence.  The  ghosted  deephole  traces  are  not  easy  to  interpret. 
Since  deghosting  the  traces  involves  amplification  at  frequencies  where  the 
signal  has  been  attenuated  by  interference  of  incident  and  surface  reflected 
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waves,  single  channel  deghosting  will  amplify  a  considerable  amount  of  noise. 
This  is  especially  evident  for  horizontals  where  excessive  ringing  is  the  result 
of  small  deviations  of  waveforms  from  the  theoretically  assumed  one  (which 
can  also  be  considered  noise) . 

If,  as  in  our  case,  several  traces  are  available,  multichannel  deghosting  seems 
to  be  a  more  attractive  approach.  It  can  avoid  increasing  the  noise  at  certain 
frequencies  by  taking  the  contribution  at  that  frequency  from  traces  where  the 
signal-to-noise  (S/N)  ratio  is  higher. 

Figure  13  shows  several  multichannel  deghosted  traces  for  a  Nicaragua  event.1 
Multichannel  deghosting  filters  proposed  by  the  Project  Officer  also  satisfy 
the  above  requirements. 

The  top  trace  shows  the  desired  output  which  is  the  actual  near  surface 
vertical  component.  The  next  two  traces  are  obtained  by  multichannel  maximum 
likelihood  filters.  The  bottom  trace  was  obtained  by  Wiener  filtering  of  the 
two  deepest  radial  components.  For  the  design  of  this  filter,  the  incident 
P  wave  was  considered  signal  and  the  reflected  SV  and  P  as  noise.  The  filter 
was  designed  to  reject  SV  and  reflected  P.  The  two  noise  components  were 
assumed  to  be  independent. 

For  maximum  likelihood  filters,  there  still  seems  to  be  some  difficulty  with 
the  inverse  filters  on  radials,  although  the  filter  no  longer  rings.  The 
Wiener  filters  are  not  optimized  with  respect  to  deephole  noise,  while  the 
maximum  likelihood  filters  are. 


3.2  ANALYSIS  OF  SEISMIC  NOISE 

In  order  to  study  the  general  character  of  noise  at  Grapevine,  several  noise 
samples  have  been  analyzed.  Power  spectra  and  cross  spectra  were  computed 
for  various  combinations  of  traces.  For  analysis  of  noise  between  10  and 
0.5  second  periods,  3000  point  samples  were  used;  128  lags  were  considered 
for  the  correlation  functions.  The  length  of  these  samples  is  2.5  minutes. 

The  horizontal  components  for  these  analyses  were  not  oriented  to  the  customary 
north  and  east,  directions,  the  component  HI  is  oriented  south-southeast, 

H2  west-southwest,  but  this  does  not  affect  the  analysis. 

Figure  14  shows  the  noise  power  spectra  computed  for  all  of  the  triaxial  components. 
Inspection  of  these  plots  shows  that  the  horizontal  noise  spectra  are  quite 
different  from  the  vertical  spectra,  which  is  to  be  expected.  The  spectra  on 
the  two  perpendicular  horizontals  are  very  similar.  The  surface  is  very 
noisy  on  verticals  and  horizontals  as  compared  to  instruments  located  at  deeper 
levels.  Surface  noise  power  around  0.5  sec  is  as  much  as  40  dB  above  the  noise 
power  level  at  the  maximum  depth  of  the  array.  The  high  surface  noise  level 
is  partially  due  to  local  geology,  near  surface  low  velocity  layer  which 


Recorded  on  2  October  1967  (location  11.7  N  latitude,  86.8  W  longitude, 
origin  time  15  h  59  m  43.4  s,  h  =  39  km  USGS) 
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Figure  14.  Triax  array  noise  power  spectra 


facilitates  the  channeling  of  seismic  energy  close  to  the  surface,  and 
partially  to  the  proximity  of  Dallas  which  is  a  source  of  cultural  noise. 


This  noise  sample  showed  very  little  directionality;  the  coherence  between 
perpendicular  components  of  motion  was  found  to  be  very  low. 

The  coherences  between  parallel  components  of  motion  at  various  depths  were 
very  low  below  1  second  period  while  between  1  and  10  second  period  higher 
coherence  values  were  found. 

A  set  of  multiple  coherences  were  computed  for  three  sets  of  outputs;  all  of 
the  HI,  H2,  and  vertical  components.  Multiple  coherence  is  defined  by  Bendat 
and  Piersol,  1966,  as 

C.2(f)  =  1  -  [Gi(f)Gi(f)] 

where  G^(f)  is  the  i7th  diagonal  element  of  the  spectral  matrix  G  of  inputs 
to  be  analyzed  and  G1 (f)  is  the  corresponding  element  of  the  matrix  G“l. 

If  n  inputs  are  being  analyzed,  n  multiple  coherences  are  obtained.  Ci  denotes 
the  relative  amount  of  power  in  the  i-th  input  trace  which  can  be  obtained 
from  the  rest  of  the  inputs  by  linear  multichannel  filtering.  Its  value  is 
unity  if  the  i-th  input  can  be  predicted  entirely  at  the  frequency  specified 
for  the  rest  of  the  inputs,  and  is  zero  if  it  cannot  be  obtained  at  all. 

Multiple  coherences  are  indicative  of  the  number  of  independent  noise  processes 
present  if  all  the  inputs  contain  these  processes,  as  in  the  case  of  a  vertical 
array  where  the  motion  at  each  level  can  be  considered  a  superposition  of 
independent  surface  wave  modes.  Theoretically,  if  the  number  of  processes  is 
less  than  the  number  of  inputs,  the  multiple  coherence  should  be  unity. 

"Independent  process"  in  the  context  of  a  vertical  array  may  mean  a  mode  with 
amplitude  depth  relationship  which  is  significantly  different  from  the  rest 
of  the  components  and  is  not  correlated  with  them.  Several  modes  with  roughly 
the  same  depth  behavior  may  behave  as  one  single  process.  The  values  of 
multiple  coherence  depend  also  on  the  relative  amplitudes  of  independent 
components.  Multiple  coherences,  therefore,  should  be  used  with  some  caution 
in  interpreting  vertical  array  noise  data.  Low  multiple  coherences  indicate 
that  an  excessive  number  of  independent  components  (modes)  contribute  significantly 
to  the  observed  inputs . 

The  squared  multiple  coherences  computed  are  shown  in  figure  15.  For  the 
vertical  components,  the  multiple  coherences  are  high  between  periods  of 
1  to  10  seconds.  Below  the  1  second  period,  all  of  the  curves  decline  steeply 
and  remain  at  very  low  coherence  values.  The  horizontal  curves  behave  similarly, 
but  they  decrease  fast  at  a  longer  period  (1.5  second).  The  higher  multiple 
coherences  for  verticals  can  be  explained  by  the  fact  that  only  Rayleigh  modes 
contribute  to  the  vertical  motion  (a  small  number  of  independent  components 
is  present).  Low  values  (below  1  second  period)  mean  that  linear  filtering 
schemes  are  not  very  effective  at  these  periods  and  the  best  array  processing 
possible  is  some  kind  of  weighted  sum.  These  low  values  indicate  that 
decomposition  of  noise  into  different  modes  is  not  possible  at  these  periods. 

The  only  mode  which  can  positively  be  identified  below  1  second  is  the 
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Figure  15.  Multiple  coherences  of  computer  transformed  horizontal 
and  vertical  channels  of  the  four  short-period  triax 
seismographs  at  GVTX 


fundamental,  causing  high  amplitudes  on  the  near  surface  instruments. 

In  general,  the  lowest  multiple  coherences  are  obtained  if  the  surface  instru¬ 
ments  are  assumed  as  output.  This  is  caused  by  the  fundamental  mode  surface 
waves  which  decrease  fast  with  depth  and  do  not  correlate  with  noise  at  the 
deepest  levels. 

Comparing  the  coherences  below  1  second  prior  ;o  theoretical  coherences  for 
P  and  S  waves  incident  uniformly  within  a  cone  of  45  degrees  aperture  at  the 
bottom  of  the  array  did  not  reveal  any  similarity.  A  previous  study  revealed 
such  a  similarity  at  Apache,  Oklahoma,  indicating  the  presence  of  mantle  P 
waves  in  the  noise  (Douze,  1967).  At  Grapevine,  this  part  of  the  short-period 
noise  seems  to  be  more  complicated. 

To  obtain  better  estimates  of  noise  spectra  at  periods  higher  than  1  second, 
the  noise  has  been  bandpass  filtered  and  decimated  by  a  factor  of  five.  A 
3600  point  sample  with  time  length  of  about  15  minutes  was  taken  and  a 
set  of  spectra  was  computed. 

These  spectra  have  a  much  higher  resolution  and  they  show  the  similarity  between 
the  two  horizontal  spectra  anu  differences  from  the  vertical.  There  is  some 
directionality  in  the  noise  around  the  6  second  period  showing  that  Rayleigh 
waves  arrive  predominantly  from  the  northeasterly  direction,  p.  ssumably 
generated  by  storms  located  in  the  Newfoundland  area. 

To  obtain  some  picture  about  the  composition  of  the  noise,  the  amplitude 
behavior  of  noise  on  the  horizontal  and  vertical  components  has  been  compared 
to  thsoretical  Rayleigh  and  Love  wave  amplitude -depth  curves. 

The  velocity  model  used  to  compute  the  amplitude-depth  curves  for  Grapevine, 

Texas,  is  given  in  table  2.  For  the  long-period  (2-10  sec)  part  of  the  com¬ 
putation,  a  rough  crustal  model  was  used  which  is  a  scaled  down  version  of 
the  model  of  Trvggvason  and  Qualls  (1966).  For  computing  the  shorter  period 
curves,  the  layered  model  is  terminated  by  a  halfspace.  The  added  crustal 
section  is  not  expected  to  influence  the  results  decidedly. 

There  is  some  uncertainty  in  the  amplitude  depth  curves  since  the  geology 
below  the  array  is  not.  well  known.  The  igneous  basement  rock  is  estimated 
to  be  about  1  km  below  the  deepest  seismometer.  If  no  drastic  changes  occur 
in  velocity,  a  few  kilometers  below  the  array,  the  curves  are  a  reasonable 
approximation.  The  densities  were  assumed  appropriate,  considering  the  lithology. 
The  S  velocities  are  computed  from  the  P  velocities  by  assuming  an  average 
Poisson  ratio  a  =  0.3  and  multiplying  the  P  velocities  by  the  appropriate 
constant  factor 


/  1  -  2a 

V  2(1  -  a)  . 

The  P  velocities  were  obtained  from  a  velocity  log. 

In  an  attempt  to  interpret  the  noise  data  in  terms  of  the  first  few  Rayleigh 
and  Love  modes,  several  deephole  to  surface  amplitude  ratios  have  been  plotted 
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Table  2.  Velocity  model  for  Grapevine,  Texas,  with  surface 
taken  as  zero  reference  and' the  first  layer  137  m  thick,  etc 


Layer 

Thickness 

(km) 

P  velocity 
km/sec 

Density 

1 

0.137 

2.43 

2.0 

2 

0.061 

2.75 

2.0 

3 

0.081 

3.50 

2.3 

4 

0.053 

2.88 

2.0 

5 

0.137 

3.50 

2.3 

6 

0 , 173 

3.04 

2.0 

7 

0.091 

3.35 

2.0 

8 

0.137 

3.50 

2.3 

9 

0.259 

3.65 

2.3 

10 

0.563 
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and  compared  to  the  corresponding  theoretical  amplitude  depth  curves.  No 
attempt  was  made  to  indicate  heavy  modes  or  body  waves  in  the  interpretation, 
although  they  may  be  present. 

The  observed  ratios  of  the  deephole  vertical  to  surface  vertical  and  the 
average  deephole  vertical  to  surface  vertical  were  superimposed  on  the  plots 
of  Rayleigh  amplitude  depth  curves.  The  ratios  of  the  averaged  deephole 
horizontals  to  the  averaged  surface  horizontals  were  superimposed  on  the  Love 
wave  amplitude-depth  curves. 

The  averaged  amplitude  ratios  for  horizontals  were  obtained  by  averaging  the 
horizontal  power  spectra,  taking  the  appropriate  ratio,  and  finally,  taking 
the  square  root  of  the  result.  This  method  approximates  the  effect  of  complete 
isotropy  in  the  noise.  In  plotting  the  ratios,  the  relative  phase  angles 
(0  degrees  or  180  degrees)  were  also  considered.  Phase  angles  were  only 
considered  if  the  coherences  between  the  respective  components  were  relatively 
high  (0.3  and  up). 

The  theoretical  horizontal  amplitudes  on  the  plot  were  reduced  by  a  factor 
l/\[T  relative  to  the  ratios  obtained  from  the  Rayleigh  wave  program  to 
simulate  the  effect  of  isotropy  of  noise. 

Plotting  the  amplitudes  with  positive  or  negative  signs  is  not  entirely 
correct  since  both  in  phase  and  180  degrees  out  of  phase  components  may  be 
prjsent.  The  plots  are,  however,  suitable  for  preliminary  analysis.  Figure 
16  shows  the  plot  obtained  for  the  6  second  period.  For  this  period,  there 
is  relatively  little  difference  in  the  Rayleigh  amplitude  depth  relationship 
on  the  vertical  components;;  the  horizontal  ratios,  however,  differ  considerably. 
The  amplitude  depth  behavior  of  the  first  three  Love  modes  is  practically 
identical. 

The  observed  vertical  ratios  are  the  closest  to  the  curve  for  fundamental 
Rayleigh  modes.  The  horizontal  V  surface  vertical  ratios  are  however,  too 
high  for  isotropic  Rayleigh  :nly.  The  observations  can  be  explained 

by  the  presence  of  fundamental  .,  leigh  waves  mixed  with  indifferentiated 
Love  modes.  Accuracy  of  the  data  does  not  permit  us  to  say  whether  any 
higher  mode  Rayleigh  waves  are  present  in  the  noise. 

Figure  17  depicts  a  similar  plot  for  the  noise  at  the  4  second  period.  The 
picture  resembles  the  previous  illustration,  but  the  vertical  ratios  are 
closer  to  the  first  higher  mode  curve.  The  fundamental  Love  mode  curve 
begins  to  separate  from  the  next  three  higher  modes.  The  observations  can 
be  explained  here  by  the  presence  of  first  higher  mode  Rayleigh  waves  mixed 
with  (preferably  higher)  Love  modes. 

Figure  18  shows  the  amplitude  ratio  plot  for  noise  at  the  period  of  3  seconds. 
The  vertical  ratios  on  this  figure  follow  the  first  higher  Rayleigh  mode 
fairly  closely.  The  horizontal  ratios  behave  quite  differently;  the  average 
phase  between  the  deepest  and  the  rest  of  the  seismometers  is  approximately 
180  degrees  and  the  horizontal  amplitudes  are  again  too  large.  The  curve 
for  the  horizontal  ratios  is  not  smooth,  an  error  resulting  probably  from 
the  fact  that  both  in  phase  and  180  degrees  out  of  phase  components  are 
present  between  the  bottom  of  the  hole  and  the  near  surface  seismometer. 
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Figure  16.  Observed  and  theoretical  amplitude  ratios 
for  noise  at  GVTX,  T  «  6,0  sec 
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Figure  17.  Observed  and  theoretical  amplitude  ratios 
for  noise  at  GVTX,  T  =  4.0  sec 
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Figure  18.  Observed  and  theoretical  amplitude  ratios 
for  noise  at  GVTX,  T  =  3.0  sec 


-29- 


Q  3707 


G  4170 


TR  68-24 


The  dominant  modes  on  horizontals  are  probably  the  higher  Love  modes  with  a 
smooth  amplitude  depth  curve  which  can  be  obtained  by  shifting  the  points 
slightly  as  indicated  on  the  figure. 

Figure  19  shows  a  similar  situation  for  a  2  second  period.  Some  Rayleigh 
modes  at  this  period  have  a  large  horizontal  amplitude  and  may  bo  contributing 
to  the  noise. 

For  the  noise  at  1  second  periods,  the  only  modes  which  can  be  identified  with 
some  confidence  are  the  fundamental  Love  and  Rayleigh  and  this  is  illustrated 
in  figure  20.  Coherences  between  levels  are  low  and  the  phase  angle  is 
uncertain. 


3.3  DEVELOPMENT  OF  PROCESSORS 

Work  has  been  done  on  two  types  of  processors  for  the  short -period  triaxial 
array.  Each  is  explained  in  some  detail  below. 

a.  A  REMODE  type  processor  was  developed  that  tests  the  character  of 
spectral  ratios  inside  a  time  window  on  all  deephole  traces.  This  processor 
chooses  one  trace  as  a  reference  and  spectral  ratios  are  computed  between  the 
references  and  each  of  the  other  deephole  traces.  The  number  obtained 
constitutes  a  complex  vector  “K.  The  normalized  complex  "dot"  product  of  vector 
A  and  the  theoretical  amplitude  ratio  for  incident  waves,  vector"^  results 
in  factors: 

«(f) 


These  factors  can  be  used  as  a  measure  of  the  similarity  of  the  motion  to  P 
motion  at  each  frequency.  The  factors  are  used  to  weight  the  Fourier 
components  of  the  reference  trace.  The  weighted  spectra  are  transformed 
point  by  point  back  to  the  time  domain  and  added  together  continuously  as  the 
window  moves  along  the  traces. 

The  method  was  applied  to  the  verticals  1,  2,  and  3  with  radial  1  as  reference; 
however,  other  combinations  of  vertical  and  horizontal  traces  can  also  be  used. 
The  output  traces  show  som*  signal-to-noise  improvement,  especially  at  the 
shorter  periods  (figure  21)  but  it  is  minor.  It  seen.3  that  at  the  longer 
periods,  the  signal  and  noise  is  very  similar  (the  deephole  array  moves  as  a 
single  unit)  and  the  parameters  are  not  sensitive  enough  for  improving  the 
signal-to-noise  ratio.  Redefinition  of  the  parameter  a  can  probably  improve 
the  method. 

Due  to  the  contract  schedule  and  the  slowness  of  this  processing  scheme,  the 
work  on  processors  based  on  similar  principles  was  not  investigated  further. 

b.  A  maximum  likelihood-REMODE  processor  is  the  second  data  processor 
developed.  Since  similarly  oriented  components  show  fairly  high  coherency 
between  various  levels,  it  is  conceivable  that  linear  processors  will  offer 
some  improvement  over  straight  summing.  Even  if  the  coherence  is  low,  linear 
processors  can  give  an  optimum  way  of  summing  (optimum  weighting)  and 
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Figure  19.  Observed  and  theoretical  amplitude  ratios 
for  noise  at  GVTX,  T  ■  2.0  sec 
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Figure  20.  Observed  and  theoretical  amplitude  ratios 
for  noise  at  GVTX,  T  =  1.0  sec 
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deghosting  the  traces .  A  REMODE  type  processor  can  then  be  used  to  accentuate 
linear  motion  on  the  deghosted,  optimum  summed  traces  that  were  derived  from 
the  group  of  radials  and  the  group  of  verticals.  Deghosting  is  necessary  to 
make  the  motion  linear  and  the  REMODE  process  is  expected  to  enhance  rectilinear 
motion  since  the  noise  between  the  horizontal  and  vertical  traces  was  found 
to  be  incoherent.  The  maximum  likelihood  filters  applied  here  are  influenced 
by  the  fact  that  the  signal  shape  is  not  identical  on  all  traces. 

The  problem  can  be  formulated  in  the  frequency  domain  as  follows:  It  is 
desired  to  minimize  the  noise  output  power  P  of  the  filters  F, 

where  F  is  a.  multichannel  filter  (complex  row  vector), 

P  =  FPnF*'  (1) 

and  PN  is  the  noise  power  spectral  matrix.  The  asterisk  denotes  complex 
conjugates  and  the  prime  denotes  transpose. 

The  filters  are  subject  to  the  condition, 

FG*  =  1  (2) 

where  G  is  a  complex  vector  denoting  the  response  of  "ghosting"  filters.  The 
condition  is  equivalent  to  stating  that  we  want  to  have  a  deghosted  waveform 
passed  by  the  filters. 

Transposing  condition  (2)  and  taking  the  complex  conjugate,  one  obtains 

G*F* '  =  1 

We  want  to  minimize  the  quantity 

FPnF*'  +  AG*F*'  where  A  is  a  Lagrange  multiplier. 

Differentiating  with  respect  to  F*'  leads  to  the  equation 

FPn  +A  G*  =  0  A  =  A/2 

multiplying  this  on  the  right  with  P  ^  and  rearranging  one  obtains 

F  =  -  AG*Pj^ 

substituting  into  condition  (2) 

FG'  =  =  1 

N 

rearranging,  one  obtains 

A  -  -1 

G*P” IfM 
N 
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substituting,  one  arrives  at  the  final  result 


G*Pjj1 

G*P^1Gi  (3) 

Equation  (3)  reduces  to  the  form  given  by  Capon  and  Greenfield  (1965)  if  the 
components  of  G  are  unity  (signal  waveform  is  identical  on  all  traces).  This 
formula  was  used  to  determine  the  spectral  response  of  the  maximum  likelihood 
filters. 

The  multichannel  filters  were  applied  to  two  groups  of  traces,  one  consisting 
of  all  radial,  the  other  of  all  vertical  components.  For  the  verticals,  the 
ghosting  filters  were  assumed  here  to  be  the  simple  (no  phase  shift)  amplitude 
response  of  form  F(u>)  *  0.5  cos  wt  where  t  is  the  uphole  time  delay.  For 
designing  the  filters,  the  quad  spectra  were  neglected  in  the  noise  spectra 
matrix.  This  is  justified  in  general  since  these  elements  were  found  to  be 
small  compared  to  the  elements  in  the  real  part  of  the  spectral  matrix. 

These  assumptions  led  to  simple,  no  phase  shift  digital  filters  which  perform 
well'  for  a  fairly  wide  range  of  the  angle  of  incidence.  The  filtering  has 
been  performed  in  the  time  domain. 

The  precursors  resulting  from  the  use  of  two  sided  filters  are  found  to  be 
small.  As  an  experiment,  all  cospectra  were  neglected  for  the  design  of 
some  filters .  This  did  not  result  in  an  appreciable  deterioration  of  the 
performance  of  the  filters  indicating  that  the  noise  is  essentially  uncorrelated 
between  the  traces  which  are  the  most  heavily  weighted  by  the  process  in  the 
frequency  range  of  interest.  The  amplitude  responses  of  the  filters  for 
vertical  components  are  shown  on  figure  22.  The  figure  shows  that  for  long 
periods,  all  traces  contribute  to  the  output;  at  short  periods,  the  contributions 
from  the  near  surface  traces  are  negligible  due  to  the  relatively  high  near 
surface  noise.  The  surface  vertical  could  be  easily  eliminated  from  the 
processing  scheme. 

The  ghosting  filter  responses  for  horizontal  components  were  found  by  the 
use  of  Haskell's  method.  The  upward  curvature  of  the  raypath  at  Grapevine  is 
due  to  the  near  surface  velocity  decrease.  The  horizontal  signal-to-noise 
ratio  is  considerably  lower  on  the  two  instruments  nearest  to  the  surface. 

The  filters  resulting  from  the  computation  have  a  phase  shift  and  the  time 
domain  representatives  are  no  longer  symmetric.  The  amplitude  responses  of 
the  filters  for  horizontals  are  shown  on  figure  23. 

Figure  24  shows  the  result  of  the  application  of  these  filters  to  an  event 
in  Nicaragua.  The  noise  sample  used  in  the  design  of  filters  was  amplified 
and  added  to  the  event  to  create  a  synthetic  noisy  record.  The  groups  of 
vertical  and  radial  components  have  been  processed  using  the  maximum  likelihood 
filters.  The  resulting  two  traces  are  shown  in  figure  24.  Also  shown  in  this 
figure  are  the  results  of  processing  the  Nicaragua  event  with  noise  added  to 
the  input  artificially  and  the  final  step  is  a  simple  motion  product  processor 
of  the  Shimshoni-Smith  type. 
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Figure  22.  Multichannel  deghosting  filter  amplitude 
response  —  vertical  components 
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Figure  24.  Maximum  livelihood  -  motion  product 
processing  Nicaragua  event 
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The  gain  factor,  (essentially  a  correlation  coefficient) 


f 


t2 

J"  v(t)«r(t)»dt 
*1 


(positive  only) , 


operates  on  the  maximum- likelihood  filtered  vertical  output  v(t)  where, 
r(t)  is  the  maximum  likelihood  radial  output,  and  t2  -  tj  is  2  seconds  (length 
of  the  moving  window).  Negative  values  of  f(t')  are  put  equal  to  zero. 


The  output  shows  that  the  noise  is  essentially  uncorrelated  between  the 
vertical  and  radial  outputs  (as  was  expected)  and  is  eliminated  quite  well. 

A  second  arrival  (likely  pP)  is  considerably  enhanced.  These  results  look 
very  promising,  but  because  of  the  poor  deghosting  achieved  on  the  horizontals, 
it  will  be  desirable  to  process  more  events  to  better  evaluate  the  method. 

The  process  is  quite  fast  and  is  probably  suitable  for  online  processing. 


4.  PROCESSING  OF  VERTICAL  ARRAY  DATA 


In  this  section  of  the  report,  a  brief  discussion  is  presented  of  all  the 
filtering  techniques  that  have  been  used  with  the  vertical  arrays  of 
vertical  seismographs.  The  triax  seismograph  discussed  throughout  this 
report  may  be  treated  as  a  vertical  seismograph  after  summation  of  the  three 
module  outputs  of  a  particular  instrument.  Therefore,  the  12-element 
triax  array  may  be  summed  to  provide  4  vertical  seismographs  operating  at 
four  depths  in  the  deephole.  In  order  to  obtain  a  coherent  presentation, 
some  portions  include  work  done  under  previous  contracts,  e.g.,  the 
section  on  inverse  filters. 

Most  of  the  time  was  employed  in  trying  techniques  based  on  least-mean- 
square  filtering.  In  general,  it  can  be  concluded  that  because  of  the  lack 
of  stationarity  of  the  noise,  these  filters  degrade  too  rapidly  to  be  of 
practical  interest  for  online  processing. 

When  the  inadequacy  of  optimum  filtering  became  apparent,  several  non-optimum 
filtering  techniques  such  as  velocity  filtering  were  tried.  In  general, 
these  non-optimum  techniques  were  as  effective  as  the  optimum  filters. 
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4.1  INVERSE  FILTERS 


One  of  the  main  differences  between  array  processing  of  deephole  and  surface 
arrays  is  that  in  surface  arrays,  the  signal  recorded  by  each  seismometer  is 
usually  assumed  to  be  identical  with  only  a  time  delay  introduced  by  the 
apparent  surface  velocity.  The  signal  recorded  at  different  depth  changes 
from  that  recorded  at  the  surface  because  of  the  interfering  reflections  from 
the  free  surface  and  velocity  interfaces.  Two  general  approaches  can  be  employed: 
(a)  the  use  of  processes  such  as  optimum  Wiener  filtering  with  assumptions 
about  the  spectrum  of  signals,  in  which  case  the  signals  do  not  need  to  be 
changed,  (Burg,  1964)  or  (b I  processes  such  as  maximum  likelihood  filters  or 
the  techniques  discussed  later  in  this  report  in  which  case  the  signals  must 
be  inverse  filtered.  We  shall  call  the  filters  that  change  a  signal  recorded 
at  one  depth  to  the  signal  at  some  other  depth  inverse  filters.  The  simplest 
approach  to  inverse  filtering  is  to  filter  the  deephole  signals  in  such  a  way 
as  to  make  them  look  like  the  surface  (deghosting  filters);  however,  it  is 
also  possible  to  change  the  signal  from  any  depth  into  the  signal  that  would 
be  recorded  at  another  depth  (ghosting  filters).  Throughout  the  present  work, 
the  filters  are  designed  so  that  all  signals  are  normalized  to  the  surface 
amplitude.  The  filters  are  designed  to  operate  on  signals  from  teleseismic 
distances  between  60  and  90  deg. 

.  The  filters  can  be  expressed  as  rational  functions  of  z,  where  the  z  transform 
(Jury,  1964)  method  is  employed. 


F(z) 


Mil  »o  ♦  *  »2Z~2  * 

»(*)  ’  i0  *  bj.-l  *  b22-2  . 


+  an2-n 
+  bmz"m 


By  rewriting  the  formula  as  B(z).  F(z)  *  A(z),  it  is  apparent  that  F(z)  is  the 
filter  that  changes  the  time  series  B(z)  into  the  time  series  A(z).  In  terms 
of  inverse  filter!:-'^  as  applied  to  deephole  signals,  the  filter  should  change 
the  signal  recorded  at  one  depth  to  that  recorded  at  another  depth.  The  time 
series  A (z)  and  B(z)  are  the  response  that  would  be  recorded  at  the  different 
depths  if  the  incoming  signal  were  a  unit  impulse;  that  is,  the  impulse 
response  of  the  layered  half  space  taken  as  a  model  v."  the  section  penetrated 

by  the  hole. 

F(z)  can  be  expanded  by  long  division  into  a  simple  polynomial: 


■  fQ  +  fiz  ^  +  f2z  ^  +  .  .  .  +  fi.z  k  +  ,  ,  ,  . 

B(z)  0  1  * 

In  general,  the  expansion  of  F(z)  will  result  in  an  infinite  series.  In  the 
filters  used  in  this  work,  the  deephole  inverse  filters  are  stable  and  the 
coefficients  converge  to  zero.  They  converge  slowly  because  of  the  large 
reflection  coefficient  at  the  free  surface.  The  series  actually  used  was 
truncated  to  a  length  of  125  points  and  a  cosine  taper  was  applied  to  smooth  the 
results.  Instead  of  obtaining  an  inverse  filter  by  long  division  and  truncation, 
W6  also  designed  a  one-channel  optimum  Wiener  filter  on  the  same  length.  The 
difference  in  the  results  was  negligible;  therefore,  the  long  division  process 
was  employed. 
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In  order  to  obtain  shorter  filters  that  would  function  more  effectively  in 
online  processing,  it  was  decided  to  use  recursive  filters  (Shanks,  1967). 
Consider  the  following  filter,  where  X(z)  is  the  input  and  Y(z)  is  the  output: 

aQ  +  a^'1  +  a2z‘2  +....+  anz'n 

Y(z)  =  - - - r - — 

1  +  b^'1  +  b2z'2  +....+  V 

where  bg  can  always  be  made  equal  to  1  (unless  it  is  zero)  by  dividing  numerator 
and  denominator  by  bg.  By  rearranging  the  formula  can  be  changed  into: 

Y(z)  =  (a0  +  ajZ-1  +  .  .  .  .  ^z'11)  X(z) 

-z'1  Y (z)  (b  x  +  b2z"1  +  .  .  bmz-m+1) 

The  output  Y(z)  is  thus  equal  to  the  input  multiplied  by  A(z)  minus  the  output 
delayed  one  sample  interval  and  multiplied  by  B(z).  In  processing  where  the 
deephole  signals  are  usually  filtered  to  approximate  the  surface  signals,  it 
is  possible  to  go  one  step  further  and  design  a  filter  that  uses  only  the 
present  input  and  past  outputs.  This  approach  is  possible  in  the  present  case 
because  A(z),  the  impulse  response  of  the  medium  at  the  surface,  consists  of 
one  large  impulse  and  some  very  minor  reverberations.  The  incoming  pulse  as 
recorded  by  the  deepest  seismometer  is  taken  as  the  reference  in  this  example. 
The  desired  filter  is: 


C(z)  can  be  obtained  by  synthetic  division  of  B(z)  by  A(z).  Because  of  the 
nature  of  A(z),  the  result  is  very  nearly  B(z)  and  the  length  of  the  filter 
has  not  increased  because  the  values  at  lags  greater  than  m  [the  length  of 
B (z ) ]  are  extremely  small  and  can  be  neglected.  Thus,  the  feedback  filter 
presently  being  used  is: 

Y (z)  -  X(z)  -  z~1Y (z)  (b2  +  b2z_1  +....+  bmz“m+1). 

It  must  be  noted  that  the  recursive  filter  is  an  exact  inverse  filter  and  not 
an  approximation  like  the  tapered  inverse  filter  discussed  previously.  Figure 
25  b,  c,  and  d  shows  the  same  signal  inverse  filtered  by  the  ghosting  and 
deghosting  techniques  with  the  long  inverse  filters  discussed  first. 

Because  of  recent  innovations  in  the  methods  of  the  Fourier  transform,  it  is 
also  feasible  to  do  the  inverse  filtering  in  the  frequency  domain.  Cooley 
and  Tukey  (1965)  have  published  an  algorithm  which  allows  very  fast  computation 
of  the  Fourier  transform. 

The  Fourier  transform  of  the  desired  inverse  filter  is  calculated  by  transforming 
the  untapered  inverse  filter  using  a  large  number  of  values  obtained  by 
synthetic  division.  This  approach  gives  a  third  method  of  inverse  filtering 
online;  but,  it  is  not  expected  to  be  more  efficient  for  online  processing  than 
the  recursive  filter.  However,  if  more  than  one  filtering  operation  has  to  be 
done,  such  as  bandpass  and  inverse  filtering  or  optimum  and  inverse  filtering, 
they  can  generally  be  done  faster  by  operating  in  the  frequency  domain.  The 
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transform  method  essentially  computes  a  Fourier  series;  therefore,  some  care 
must  be  employed  in  using  the  method  for  filtering.  Convolutions  as  usually 
computed,  assume  the  time  series  to  be  zero  outside  the  region  of  interest. 

The  new  method  assumes  the  input  to  be  periodic  with  a  period  equal  to  the 
length  of  the  record.  When  multiplying  in  the  frequency  domain,  the  periodicity 
implies  that  elements  shifted  past  the  end  of  the  time  series  appear  again  at 
the  beginning  (personal  communication,  McCowan) .  This  problem  can  be  avoided 
by  adding  a  number  of  zeros  to  the  end,  equal  to  the  number  of  filter  points 
in  the  time  domain.  When  this  is  done,  no  end  effects  are  present  and  frequency 
domain  filtering  will  result  in  a  smooth  time  series  with  no  discontinuities 
at  the  points  where  the  time  series  divided  and  transformed  into  the  frequency 
domain  and  back.  It  is,  of  course,  possible  to  compute  the  inverse  filters 
in  the  frequency  domain.  However,  the  computer  program  used  allows  specifications 
of  the  filter  in  either  the  time  domain  or  frequency  domain  and  adds  the  correct 
number  of  zeros  or  the  number  of  zeros  needed  if  the  filter  is  defined  in  the 
frequency  domain.  In  the  frequency  domain,  the  number  of  frequency  points  must 
be  equal  for  filter  and  data;  thus,  a  large  number  of  zeros  are  usually  added 
to  the  time  domain  filter  to  make  the  length  of  filter  and  data  equal.  Inverse 
filtering  in  the  frequency  domain  resulted  in  a  deghosted  signal  identical  to 
that  shown  in  figure  25  b ,  c,  and  d. 


4.2  OPTIMUM  FILTERING 

Two  methods  of  optimum  filtering  have  been  tried.  In  the  first,  the  signals 
are  equalized  at  all  depths  by  the  use  of  inverse  filters  and  the  signals 
were  eliminated  by  subtraction;  optimum  filters  operating  only  on  the  noise 
were  then  computed.  The  second  approach  used  consists  of  standard  optimum 
least-mean-square  filters  in  which  assumptions  about  the  signal  must  be  made. 
Discussion  of  optimum  filtering  in  this  report  will  be  confined  to  Wiener 
least-mean-square  tecnnique  discussed  by  Burg  (1964) . 

4.2.1  Deghosting  Optimum  Filters 

The  processes  used  are  shown  in  figure  26  „  These  processes  are  logical 
extensions  of  a  technique  discussed  briefly  by  Bendat  (1958).  It  must  be 
noted  that  the  inverse  filters  must  be  quite  exact  for  the  process  to  operate 
correctly.  It  was  expected  that  for  large  signals,  small  errors  in  the  inverse 
filters  would  result  in  the  subtraction  shown  in  figure  26  not  being  exact  and 
errors  would  be  introduced.  This  phenomenon  is  unavoidable  when  a  range  of 
distances  and  angles  of  incidence  are  covered  by  a  single  set  of  inverse 
filters.  However,  large  signals  with  a  good  signal-to-noise  ratio  are  of  little 
interest  in  array  processing  because  all  the  necessary  parameters  can  be 
measured  from  a  single  seismograph.  As  will  be  shown  later  in  the  discussion, 
the  inverse  filters  are  sufficiently  exact  to  allow  the  process  to  operate 
correctly. 

These  processing  schemes  are  based  on  the  concept  that  elimination  of  the 
assumptions  made  about  the  signal  in  conventional  optimum  filtering  would 
result  in  a  more  effective  processor.  Figure  26  shows  the  filtering  technique 
used  in  the  first  attempts  to  improve  the  detection  capability  of  the  site. 

As  shown  by  the  schematic,  the  signals  are  first  equalized  by  the  inverse 
filters  discussed  in  the  previous  section.  After  inverse  filtering,  each 
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Figure  26.  Deephole  data  processing  techniques 


deephole  channel  is  subtracted  from  the  reference  channel  to  produce  a  trace 
that  consists  of  noise  only.  Although  the  surface  channel  has  been  used  as  a 
reference  channel  throughout  the  present  analysis,  any  channel  can  be  used 
as  a  reference.  In  fact,  at  sites  where  there  is  an  appreciable  increase 
o±  signal-to-noise  ratio  with  depth  (which  is  not  the  case  at  APOK),  the 
surface  must  not  be  used  as  a  reference  because  it  would  introduce  excessive 
noise  into  the  process . 

An  optimum  filter  was  then  constructed  where  the  desired  signal  is  the  noise 
at  the  surface  and  the  actual  noise  recorded  in  the  deephole  is  to  be  suppressed. 
It  must  be  noted  the  usual  assumption  that  desired  signal  (surface  noise)  and 
noise  are  uncorrelated  is  not  correct;  the  optimum  filters  are  designed  from 
actual  data  (including  correlation  of  signal  and  noise)  and  no  assumptions 
about  the  signal  are  made.  The  optimum  filter  output  is  then  subtracted  from 
the  surface  to  obtain  the  desired  signal  (surface  noise).  The  process  can 
use  either  multi-  or  single- channel  optimum  filters.  Both  the  inverse  filters 
and  optimum  filters  are  operating  correctly  as  shown  in  figure  25.  The  shape 
of  the  signals  (for  25b  through  25g)  is  identical  to  that  recorded  by  the 
surface  seismograph.  It  must  be  noted  that  this  sharp  signal  constitutes  a 
severe  test  for  the  inverse  filters  (recursive). 

Figure  27  shows  the  decrease  of  the  noise  level  obtained  by  some  of  these 
processes.  The  single -channel  filters  are  all  one-sided  filters  using  only 
past  data  because  experimental  comparisons  indicated  that  no  further  improve¬ 
ment  is  obtained  by  two-sided  filters  in  this  situation. 


The  filter  length  employed  was  71  points,  corresponding  to  3.5  sec  of  real 
time.  In  comparing  the  three  single-channel  outputs,  it  is  apparent  that  the 
noise  level  was  only  reduced  by  a  small  amount  except  at  0.5  sec.  The  reason  for 
this  behavior  is  probably  that  the  amplitude-depth  relationships  of  signal  and 
noise  are  almost  identical.  The  peak  at  3.0  sec  period  was  not  reduced  appreciably 
by  any  of  the  filters  mainly  because  the  filter  length  was  not  sufficient;  however, 
as  these  waves  are  not  at  the  periods  where  teleseisms  are  usually  detected, 
no  attempt  was  made  to  improve  the  performance . 

The  2  Hz  peak  shown  on  all  spectra  is  one  main  obstacle  to  visual  detection  at 
the  APOK  site.  Comparison  of  the  results  of  the  single-channel  cases  indicate 
that  deephole  1  (2880  m)  was  the  most  effective  and  reduced  the  power  in  the 
2  Hz  peak  by  12  dB.  The  single  channel  process  operating  on  deephole  3 
(1060  m)  did  not  reduce  the  power  while  deephole  2  (1980  m)  achieved  a  reduction 
of  5  dB.  Therefore,  a  single-channel  optimum  filter  operating  on  the  bottom 
seismograph. is  the  best  processor  for  this  technique. 


Figure  27  also  shows  the  results  of  using  a  single -channel  ghosting  technique, 
in  which  case,  all  the  recordings  were  filtered  to  simulate  deephole  1  before 
subtraction.  For  reasons  that  are  not  understood  at  present,  the  results  are 
inferior  to  the  single-channel  deghosting  filter  technique  if  the  attenuation 
of  the  spectral  peak  at  2  Hz  is  used  as  the  criteria.  At  0.5  Hz,  the  ghosting 
technique  resulted  in  a  greater  attenuation  of  the  noise.  The  process  is 
working  correctly  as  indicated  in  figure  25  in  which  deephole  1  (unfiltered) 
and  the  final  output  of  the  process  are  shown. 
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Figure  27.  Attenuation  obtained  from  optimum  techniques:  a/  single¬ 
channel  deghosting  (2880  m)  b/  single-channel  deghosting 
(1980  m)  c/  single-channel  deghosting  (1060  m)  d/  sum  of 
single-channels  e/  one-sided  multi-channel  deghosting 
f/  two-sided  multi-channel  deghosting  g/  sum  of  deghosted 


A  two-channel  Wiener  optimum  filter  set  for  GVTX  was  computed  using  the 
technique  in  which  the  signals  are  eliminated  by  subtracting  the  inverse 
filtered  traces.  DH1  was  used  as  a  reference  and  the  surface  seismograph 
was  not  used.  Only  slight  attenuation  in  the  noise  was  obtained  when  the 
filters  were  applied  to  the  noise  from  which  they  were  calculated.  When 
operating  on  noise  samples  several  days  later,  no  attenuation  of  the  noise 
was  obtained.  It  is  of  interest  to  note  that  the  signal  shape  remained 
undistorted  for  several  cycles.  This  phenomenon  may  indicate  that  signal 
generated  noise  is  not  as  much  of  a  problem  at  GVTX  as  it  is  at  APOK. 

4.2.2  Conventional  Optimum  Filters 

For  comparison  with  the  results  discussed  above,  optimum  filters  of  a  more 
conventional  type  were  also  calculated.  The  technique  consists  of  assuming 
a  model  of  the  signal.  In  the  case  discussed  here,  the  signal  is  assumed  to 
be  white  between  2.0  and  0.3  sec  period  and  a  cosine  taper  starting  at  the 
cutoff  is  zero  at  4.0  and  0.4  sec  period.  The  impulse  response  of  the 
velocity  section  in  the  hole  is  convolved  with  the  theoretical  signal  to 
obtain  the  model  of  the  signal  at  any  depth.  Signal  and  noise  are  taken  to 
be  uncorrelated.  The  correlation  matrix  then  consists  of  the  measured  noise 
added  to  the  theoretical  signal  for  the  various  depths.  The  desired  signal 
is  the  signal  as  recorded  by  the  surface  seismograph.  It  should  be  noted 
that  an  assumption  must  be  made  as  to  the  signal -to-noise  ratio.  A  signal-to- 
noise  of  2  was  assumed  in  these  calculations. 

Using  the  deghosting  techniques  discussed  previously,  it  was  found  that  a 
single-channel  filter  at  APOK  was  as  effective  as  a  multichannel  filter. 
Therefore,  it  was  decided  to  try  both  2-channel  and  4-channel  optimum  filters. 
Figure  27  shows  the  noise  reduction  that  was  obtained;  if  the  2  Hz  peak  is  again 
used  as  a  criteria,  these  optimum  filters  are  as  effective  as  the  single¬ 
channel  filter  using  deghosting  technique  in  reducing  the  noise.  Examination 
of  figure  27  indicates  that  the  2-channel  filter  is  more  effective  than  the 
4-channel  filter.  Examination  of  figure  28  illustrates  that  the  filter  also 
reduces  the  signal  amplitudes  approximately  3  dB.  It  is  of  course  possible 
to  use  a  number  of  different  signal  models  and  signal-to-noise  ratios  to  try 
to  improve  on  the  results  obtained  here.  However,  it  would  appear  likely 
that  an  improvement  in  noise  reduction  would  result  in  a  further  reduction  of 
signal  amplitude. 

The  results  obtained  to  date  indicate  that  the  one-sided  single-channel  optimum 
filter  using  the  deephole  1  and  the  surface  (as  reference)  is  as  effective  a 
technique  as  any  multichannel  set  of  filters  at  APOK.  However,  because  of 
the  complexity  of  programming  the  first  technique  online,  it  is  probably 
preferable  to  use  conventional  multichannel  filters. 

A  three-channel  optimum  Wiener  filter  was  computed  using  the  three  deepest 
seismometers  at  GVTX.  The  surface  seismograph  was  not  included.  In  theory, 
inclusion  of  the  high  noise  level  at  the  surface  would  result  in  a  set  of 
filters  which  would  heavily  attenuate  the  signal  contribution  from  the  surface 
seismograph.  In  practice,  the  filter  could  not  hope  to  respond  to  the  large 
burst  of  noise  from  nearby  traffic  without  degrading  the  performance  from 
the  array.  The  signal  model  used  consists  of  a  white  signal  between  0.25 
and  2.0  second  period,  with  a  cosine  taper  at  each  end.  This  model  appears 
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to  cause  less  trouble  in  computation  of  filters  than  the  model  without  the 
tapers.  The  results  are  shown  in  figure  29,  together  with  the  ncise  level 
of  DH1,  which  has  the  best  signal-to-noise  ratio  of  the  individual  seismographs 
Only  a  small  amount  of  noise  reduction  was  obtained. 

4.2.3  Prediction  Filters 

A  set  of  optimum  prediction  filters  were  computed  for  the  APOK  vertical  seismo¬ 
graphs.  These  filters  do  not  predict  noise  in  the  future,  but  take  advantage 
of  the  fact  that  the  signal  has  already  been  received  by  the  bottom  seismograph 
while  the  other  seismographs  are  still  recording  noise.  Thus,  a  set  of  filters 
operating  on  the  noise  from  the  top  seismographs  can  be  used  to  predict  the 
noise  from  the  bottom  seismograph^  subtracting  the  two  time  series  should 
eliminate  some  of  the  noise  and  improve  the  signal-to-noise  ratio  of  the  first 
break.  Of  course,  it  must  be  noted  that  when  the  signal  arrives  at  the  other 
seismographs,  the  signal  will  be  distorted.  Because  of  the  similarity  between 
noise  and  signal  at  the  site,  the  signal  is  severely  attenuated  after  the 
first  break.  Figure  30  shows  the  noise  after  3-channel  prediction  filtering 
and  the  unfiltered  noise  at  the  surface.  In  this  case,  the  noise  reduction  is 
greater  with  increasing  numbers  of  seismographs.  As  noted  previously,  optimum 
filtering  using  other  techniques  does  not  always  result  in  further  improvements 
when  more  seismographs  are  added. 

It  was  decided  to  try  long  prediction  filters  (>200  points,  20  samples/sec) 
to  determine  if  further  improvements  are  possible.  The  total  amount  of  data 
used  was  3276  sec  of  noise  with  a  spectral  lag  of  12.8  sec.  The  results  are 
high  confidence,  high  resolution,  spectral  and  cross  spectral  estimates.  The 
multiple  coherences  were  calculated  to  determine  the  improvement  that  can 
theoretically  be  obtained.  The  multiple  coherence  is  defined  by  Bendat  (1966): 

Y-j  (£)  -  1  -  [Gi(f)  Gi(f)]  1 

where,  G1(f)  is  the  ith  diagonal  element  of  the  inverse  of  the  spectral  matrix. 
The  noise  reduction  that  can  theoretically  be  obtained  is: 

y«  -f1  -  v-xW« 

where  G^(f)  is  the  theoretically  predicted  spectrum.  , 

Figures  31  and  32  show  the  multiple  coherences  at  APOK  when  the  surface  and 
DH1  (2880  m)  are  the  outputs.  The  multiple  coherences  are  high  for  periods 
greater  than  1.0  sec  and  decrease  rapidly  towards  the  shorter  periods.  When 
the  surface  is  taken  as  the  output,  the  2  Hz  peak  shows  a  coherence  somewhat 
higher  than  the  adjacent  frequencies.  With  DH1  as  the  output,  the  coherence 
at  the  2  Hz  peak  is  no  more  than  adjacent  frequencies.  In  general,  the  plot 
shows  the  higher  frequency  noise  is  not  very  predictable  when  averaged  over 
long  samples. 

A  210  point  set  of  3-channel  prediction  filters  were  computed  in  the  time 
domain  from  the  correlation  matrix  obtained  from  the  long  sample  discussed 
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Figure  29.  Noise  spectrum  of  DH1  in  solid  line,  noise  spectrum  after 
3~channel  optimum  filtering  in  dashed  line,  GVTX 
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Figure  30.  Spectrum  of  the  noise  after  three-channel  prediction 
filter  using  DH4,  DH3,  and  DH2  to  predict  DH1  (solid 
line),  and  spectrum  of  the  noise  at  the  surface 
(dashed  line)  _5Q_ 
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Figure  32.  Multiple  coherence  of  the  four  vertical  seismographs  at 
APOK,  DH1  used  as  output 


before  Figure  33  shows  the  theoretically  predicted  noise  spectra  after 
prediction  filtering  and  the  actual  spectra  obtained  experimentally  by  applying 
the  set  of  filters  to  the  first  240  sec  of  the  data  from  which  the  filters 
were  computed.  At  no  frequency  are  the  experimental  results  as  good  as  the 
multiple  coherence  indicates  they  should  be.  At  2  Hz,  no  improvement  was  obtained. 
It  is  concluded  from  these  experiments  that  the  noise  is  not  stationary  at  2  Hz 
over  periods  on  the  order  of  an  hour.  At  1  Hz,  the  noise  is  not  completely 
stationary,  but  some  improvements  can  be  obtained  by  optimum  filtering. 

In  order  to  compute  long  prediction  filters,  a  program  has  been  written  to 
obtain  spectral  and  cross  spectral  estimates  from  large  noise  samples.  The 
program  uses  the  Fast  Fourier  Transform  and  operates  by  sectioning  the  record, 
taking  the  spectra  and  averaging  these  spectra. 

Let  X(t) ,  t  «  0,  .  .  .  .  N-l  be  a  sample  from  a  seismic  time  series.  We  take 
non- overlapping  segment  X^Ct),  t  -  0,  1,  .  .  .  L-l,  where  n(L-l)  =  N-l. 

The  Fourier  transform  of  each  segment  is  computed 

L- 1 

Fn(k)  =  i  ^  Xn(t)  w<t)e'2lTikn/L 
t=0 

where  W(t)  is  the  smoothing  window;  to  date  a  Hanning  window  has  been  used. 

Fn(k)  *  1/AF^i  *  l/2Fn  *  l/4Fn+1 

The  spectral  and  cross  spectral  of  each  segment  are: 

S  (k)  =  F.  (k)  F*(k) 
n  i  j 

where  i*j  gives  the  spectra  and  i/j  gives  the  cross  spectra. 

The  final  step  consists  of  averaging  over  all  the  spectra  to  obtain  the  final 
spectral  estimate  of  X(t), 

sew  - 1  2  Vk> 

n=0 

The  average  of  this  method  consists  of  not  being  limited,  in  the  length  of 
sample  that  can  be  analyzed,  by  the  core  of  a  small  computer.  For  example, 
we  have  calculated  the  spectral  matrix  from  1  hour  of  short-period  data  with 
256  lags  in  only  25  minutes  of  CDC  3200  computer  time.  The  program  is 
specially  useful  in  obtaining  high  resolution  and  high  confidence  multiple 
coherences . 

4.3  NON-OPTIMUM  FILTERS 

As  an  alternative  to  the  use  of  optimum  filters,  which  deteriorate  rapidly  with 
time  because  of  the  non- stationary  noise,  a  number  of  non-optimum  filtering 
techniques  (in  the  least-mean-square  sense)  were  tried.  The  methods  tried 
during  the  present  contract  were  multichannel  deghosting  filters,  velocity 
filters,  and  time  delay  and  sum. 
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Figure  33.  Experimental  results  from  prediction  filtering  (solid 
line)  and  theoretically  predicted  spectrum  after 
prediction  filtering,  (dashed  line)  APOK 
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4.3.1  Multi cnannel  Deghostir.g  Filters 


A  set  of  multichannel  filters  were  computed  from  a  formula  suggested  by  the 
Project  Officer. 


M  * 
O 

t — \ 

W  1 

a  cos  2irft 

0  o 

zx(f) 

— 

cos  2irftj 

zn(f)_ 

an  cos  2irftn 

*S(f) 


where  n  is  the  number  of  seismographs,  tn  is  the  uphole  time,  S(f)  is  the 
surface  signal  and  a  weighting  factor  obtained  from  the  noise  level  of  the 
nth  seismometer.  The  filter  applied  to  each  channel  is, 
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The  advantage  of  this  filter  is  that  the  noise  level  cannot  increase  above 
the  level  of  an  individual  seismograph.  The  filters  were  computed  in  the 
frequency  domain  and  transformed  into  the  time  domain.  Filters  obtained  in 
this  fashion  are  two  sided  and  small  precursors  appear  on  the  filtered  seismo¬ 
grams  in  front  of  the  signal.  Figure  34  shows  the  filter  result  obtained 
using  all  four  seismographs  at  GVTX  and  the  three  bottom  seismographs.  The 
slight  difference  in  amplitudes  is  caused  by  the  presence  or  absence  of  the 
surface  trace  with  the  associated  weighting  factor.  A  small  precursor  is 
visible  on  both  examples,  but  is  more  noticeable  on  the  lower  trace.  A  com¬ 
parison  with  the  surface  trace  indicates  that  the  filters  are  operating 
correctly.  The  main  advantages  of  this  method  are  the  simplicity  of  computation 
and  the  fact  that  it  does  not  allow  for  an  increase  in  the  noise  levels.  The 
main  disadvantages  are  that  long  filters  (150  points  were  used)  are  necessary 
and  that  precursors  appear  in  front  of  the  signals. 

4.3.2  Velocity  Filters 

While  body-wave  signals  travel  uphole  with  known  velocity,  that  part  of  the 
noise  which  consists  of  Rayleigh  waves  is  either  in  phase  or  180  deg  out  of 
phase  with  the  surface.  Under  these  conditions,  it  appeared  that  velocity 
(fan)  filters  might  be  effective  in  reducing  the  noise  level. 

This  approach  was  used  only  with  data  from  APOK,  where  the  signal-to-noise 
ratios  are  close  to  constant  with  depth.  The  technique  is  not  applicable  to 
GVTX  where  the  surface  seismograph  should  not  be  used  and  therefore  only  a 
small  linear  array  is  available. 

A  simple  set  of  velocity  (fan)  filters  were  calculated  for  the  APOK  deephole 
site.  The  data  were  single-channel  inverse  filtered  and  then  fan  filtered  to 
pass  only  waves  traveling  at  the  uphole  velocity  of  P  waves.  The  filters  used 
were  calculated  from  the  following  formula  (Embree,  et  al,  1963): 
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(b)  bottom  channels,  GVTX 


1 


where:  Tn  =nAt  a  given  sampled  value  of  time  t 

Ax  ^seismograph  separation 

•  til 

“distance  from  array  center  of  m  ‘  seismograph. 

Because  of  the  small  array  available,  the  filters  attenuated  both  the  signal 
and  the  noise.  Figure  35  shows  the  noise  attenuation  obtained  together  with 
a  signal  before  and  after  filtering.  Comparison  of  signal  and  noise  attenuations 
shows  that  for  periods  less  than  1.0  sec,  a  small  improvement  was  obtained. 
Improved  fan  filters  can  be  calculated  and  should  be  tried.  Finally,  it  must 
be  noted  that  with  this  very  small  array,  only  a  limited  improvement  on  the 
signal- to-noise  can  be  expected. 


4.3.3  Time -Del ay  and  Sum 

In  order  to  compare  with  the  results  of  filtering  techniques  described  in  the 
previous  sections,  it  was  decided  to  include  simple  time  delay-and-sum 
filtering.  The  major  disadvantage  of  the  method  is  that  for  vertical  arrays, 
the  resulting  signal  shape  is  not  the  one  that  is  recorded  at  the  surface. 

The  major  advantage  is  the  simplicity  of  the  operation. 


Several  different  weighting  functions  were  used  in  various  time-delay  and  sum 
operations  at  GVTX  and  are  listed  below: 


No.  1  No.  2  No.  3 


DH4 

0.149 

0.227 

0.000 

DH3 

0.373 

0.555 

0.373 

DH2 

0.612 

0.833 

0.612 

DH1 

1.000 

1.000 

1.000 

The  weighting  factors  for  example  No.  1  were  provided  by  the  Project  Officer. 
Weighting  factors  for  example  No.  2  were  calculated  from  the  mean  square 
noise  level  of  the  spectra.  The  spectra  are  corrected  to  equal  signal  level. 
Figure  36  shows  the  results  obtained  and  neither  resulted  in  an  improvement 
over  the  noise  level  at  the  bottom  of  the  hole.  Figure  36  also  shows  the 
results  from  the  weighting  factors  of  example  No.  3  where  the  noise  surface 
seismograph  was  not  included.  The  results  show  a  modest  improvement  over  the 
noise  level  at  the  bottom  of  the  deephole.  From  these  experiments,  it  appears 
that  it  is  better  not  to  include  the  surface  seismographs  at  sites  with  high 
level  surface  noise. 
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NOISE  POWER  ATTENUATION 


PERIOD  (second) 


G  3702 


Figure  35.  Noise  (from  spectra)  and  signal  attenuation  resulting  from 
fan-filtered  deephole  array,  APOK 
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It  is,  of  course,  always  possible  to  include  the  surface  by  assigning  it  a 
weighting  factor  consistent  with  the  noisj  level  when  large  bursts  of  noise 
are  present.  However,  in  this  case,  the  weighting  factor  for  the  surface  at 
Grapevine  would  be  about  0.05  and  is  therefore  essentially  negligible. 

At  APOK,  the  time  delay  and  sum  operations  were  done  without  weighting  factors 
because  the  signal-to-noise  ratio  is  almost  identical  for  all  depths.  Only 
a  small  increase  in  the  signal-to-noise  ratio  (3  dB)  was  obtained  m  this 

fashion. 
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6.  CONCLUSIONS  AND  RECOMMENDATIONS 


6.1  CONCLUSIONS 

The  short-period  triaxial  seismometer  has  proven  to  be  a  very  stable  and 
reliable  instrument,  requiring  a  minimum  of  attention  and  adjustments.  It 
is  also  much  easier  to  work  with  and  transport  than  earlier  deephole  instruments. 
Nearly  all  of  the  malfunctions  encountered  in  the  downhole  instrumentation  can 
be  attributed  to  dc  motor  failure  or  to  insulation  breakdown  in  the  7H4  cable. 

The  insulation  in  the  7H4  cable  used  on  this  contract  has  broken  down  on 
several  occasions  due  to  age  and  long  periods  of  use  in  high  temperature  deep- 
holes.  Cable  insulation  problems  were  the  cause  of  most  of  the  seismometer 
downtime. 

Considerable  time  was  lost  due  to  the  computer  hardware  problems.  As  a  result, 
online  program  development  was  delayed  during  the  beginning  of  the  contract 
year.  After  an  on-call  service  contract  was  made  with  Adage,  Inc.,  the 
computer  manufacturer,  satisfactory  operation  was  maintained.  Normal  operation 
of  the  computer  was  possible  from  the  last  week  of  August,  1967,  until  the 
Grapevine,  Texas  site  was  shut  down  the  last  week  of  May,  1968. 

None  of  the  filtering  techniques  employed  in  the  recordings  of  the  vertical 
array  have  resulted  in  improved  detection  capability  of  the  magnitude  originally 
expected.  Least-mean-square  processing  methods  are  not  acceptable  because  the 
filters  computed  from  a  given  noise  sample  degrade  rapidly  (<  1  hour)  with 
time.  Therefore,  these  techniques  are  not  usable  for  online  processing. 

It  must  be  noted  that  adaptive  optimum  filtering  could  invalidate  this  conclusion 
However,  it  appears  unlikely  that  adaptive  filtering  will  be  any  better  than 
Wiener  optimum  filters  used  on  the  sample  from  which  they  are  computed.  There¬ 
fore,  even  with  adaptive  filtering  improvements  in  signal-to-noise  of  3-4  dB 
over  the  period  range  of  interest  appears  to  be  the  best  that  can  be  expected. 

The  non-optimum  techniques  of  multichannel  deghosting  and  velocity  filtering 
are  almost  as  effective  as  optimum  filtering  and  need  only  be  computed  once  for 
each  site.  These  techniques  are  therefore  to  be  preferred  for  online  processing. 


6.2  RECOMMENDATIONS 

Dc  motors  were  the  only  components  within  the  seismometer  that  failed  and 
caused  damage  to  the  instrument  itself.  Oxidation  and  corrosion  of  the 
armature  and  brushes  appear  the  major  cause  of  motor  failure.  These  failures 
prevent  the  seismometer  masses  from  being  caged  and  consequently  the  suspensions 
are  very  often  damaged  while  removing  the  seismometers  from  the  dcsphole.  In 
order  to  correct  this  motor  problem,  it  is  recommended  that  brushless  stepping 
motors  be  tested  and  used  to  replace  the  existing  motors. 

Armored  electrical  cable  (7H4)  nsed  on  this  program  is  becoming  quite  unreliable 
because  of  normal  usage  over  several  years.  It  is  recommended  that  new  cable 
be  supplied  if  the  system  is  transported  to  another  location  for  data  acquisition 
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A  large  amount  of  valuable  equipment  is  now  in  storage  as  a  result  of  the 
termination  of  this  contract.  This  equipment  will,  in  all  likelihood,  be 
scattered  within  other  programs  if  it  is  not  put  back  in  operation  in  a 
reasonable  length  of  time.  It  is  recommended  that  this  short-period  triaxial 
system  be  used  to  replace  older  equipment  at  some  location  where  high  quality 
deephole  data  may  be  recorded  for  further  research  work. 

The  maximum  likelihood  REMODE  nonlinear  processor  discussed  in  section  3.3 
appears  quite  promising.  Because  of  the  shortage  of  time  and  funds,  this 
technique  has  not  been  fully  evaluated.  It  is  recommended  that  additional 
analysis  of  three -component  deephole  data  be  included  in  future  work. 
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APPENDIX 


REPRODUCTION  OF 

_  STATEMENT  OF  WORK  TO  BE  DONE 

AFTAC  PROJECT  AUTHORIZATION  NO.  VELA  T/7703 
AND  AMENDMENTS  _ 


APPENDIX  1  F33657-67-C-1224 

18  QtU  19bb 

SIV.'t .'A1ENT  OF  V.'OPX  TO  BE  PONE 
(AFTAC  Project:  Authorization  No.  VI/7703/S/ASD) 

1,  Usin"  the  deap-well-array  field  system  assembled  for  Project  VT/5051, 
conduct  a  fie id-opur.ic ions  program  at'  two  or  more  available  VELA-UNTFORM 
deep-well  sites  to  he  designated  by  the  AFTAC  project  officer.  The  sites 
near  Ptnedale,  Wyoming  (PIl-IY)  .  Grapevine,  Texas  (G7-TX),  and  Franklin,  West 
Virginia  (FEW/) ,  may  be  considered  for  planning  purposes.  In  the  course 
of  this  work,  accomplish  the  following  subtasks  at  each  site  visited: 

a.  Routinely  operate  <  nd  maintain  an  earily-transportable  deep-well 
vertical-array  system  ploying  4  to  6  triaxisra  and  vertical  component 
short-period  borehole  seismometers,  appropriate  surface  instrumentation, 
and  an  on-line  Jigitizing/reco-ding/computing  device. 

b.  Record  on  film  and  digital-tape  vertical  and  three-component  data 
from  the  deep-well  array.  Establish  an  interim  library  of  digitized  raw 
seismic  data  from  the  array,  including  records  of  background  noise  and 
signals  and  appropriate  identifying  logs,  suitable  for  use  in  this  and 
other  projects. 

c.  Operate  the  on-line  computing  device  to  test  in  real-time  different 
vertical-array  signal- enhancement  data-proceseing  techniques.  Record  the 
results  on  film  and  tape.  Operate  at  each  site  long  enough  to  compare  the 
different  techniques  and  determine  the  better  signal-enhancement  processing 
methods  for  each  site.  Relate  the  results  of  each  technique  to  site 
conditions,  such  as  seismic  noise  field,  geological  structure,  and  geographic 
location.  This  work  is  to  be  accomplished  in  close  coordination  with  the 
AFTAC  project  officer,  who  will  assist  in  selecting  the  techniques  to  be 
tested. 

2.  Process  the  raw  seismic  data  obtained  from  field  operations  and 
perform  detailed  analysis  and  evaluation  of  the  horizontal  components 
of  surface  and  subsurface  noise  at  the  designated  sites,  making 
comparisons  between  observations  and  theoretical  studies,  to  improve 
our  understanding  of  noise  composition.  Particular  emphasis  is  to  be 
placed  on  using  horizontal  data  to  ascertain  the  true  form  of  noise 
constituents  which  appear  from  their  vertical  components  to  be 
compressional  body-waves.  Computer  programming  and  computational  services 
will  be  accomplished  by  the  contractor;  commercial  high-speed-computer 
services  will  not  be  arranged  for  under  this  project  without  prior  AFTAC 
concurrence. 


3.  Tailor  various  borehole -array  data-processing  techniques,  for  improving 
signal-to -noise  ratios  and  detecting  low-level  signals,  for  implementation 
at  the  designated  deep-well  sites  and  demonstrate  them  in  off-line  tefcts. 
The  processes  should  be  ultimately  suitable  for  real-time  use  and  shall 
include,  but  not  necessarily  be  limited  to:  optimum  filtering,  both  with 
and  without  ghosting  and  deghosting  pre-filters;  prediction  filtering, 
using  vertical  data  exclusively  and  three-component  data;  correlation 
processing  and  other  non-linear  techniques  and  possibly,  ether  processes 
suggested  by  the  AFTAC  project  officer.  Adapt  selected  signal-enhancement 
techniques  to  on-line  real-time  operation. 


;  Uachrucnt  1 
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REPORTS 

(AFTAC  Projeci:  Authorization  No.  VT/7703) 

1.  Monthly  14a. p.2 cement  Report;  A  monthly  letter-type  management  report, 
summarizing  work  for  the  calendar  month  and  generally  not  exceeding 
3  pages  of  text,  will  be  submitted  to  AFTAC  in  15  copies  by  the  20th  day 
of  the  following  month.  The  heading  of  each  report  will  contain  the 
identification  data  lir.'ed  in  paragraph  7.  These  reports  should  not 
contain  detailed  t-_chnic.il  descriptions,  formulae,  equations,  graphs,  or 
other  technical  documentation,  instead,  they  will  present  a  narrative 
summary  of  work  performed  during  the  report  period,  including  reference 
to  the  following  topics: 

a.  Research  Program  find  Plan  ••  a  brief  statement  of  objectives  and 
plan  of  research. 

b.  Major  Accomplishments  ••  a  brief  description,  written  in  lay  terms, 
of  any  finding  or  accomplishment  considered  worthy  of  being  brought  to  the 
attention  of  management. 

c.  Problems  Encountered  -  make  reference  to  difficulties  associated 
with  personnel,  facilities,  contracts,  availability  of  literature,  funds, 
strikes,  disasters,  etc.,  which  significantly  affect  the  progress  of  the 
work  involved.  Problems  of  a  technical  nature  should  also  be  included,  but 
in  brief,  nontechnical  terms. 

d.  Fiscal  Status  -  include  an  estimate  of  contractor  funds  required 
to  complete  the  work  if  an  overrun  is  anticipated. 


e.  Government  Actions  Required  -  specify  any  AFTAC  assistance  required 
in  resolving  "Problems  Encountered." 

f.  Future  Plans  -  a  brief  statement  of  any  significant  change  which  is 
planned  in  the  course  of  work  under  way  or  on  any  new  item  which  is 
considered  to  be  of  interest  to  management.  In  addition,  the  initial 
monthly  report  should  include  an  introduction  outlining  the  background, 
objectives,  and  assignment  of  responsibility  for  the  project. 

2.  Monthly  Financial  Statement:  The  contractor  will  submit  5  copies  of 
financial  data  in  an  Alternate  Management  Summary  Report. 

3.  Quarterly  Technical  Reports:  Quarterly  letter-type  technical  reports, 
summarizing  work  for  calendar  quarters,  will  be  submitted  to  the  AFTAC 
project  officer  in  2  copies  within  15  days  after  the  close  of  each  such 
period.  Each  report  will  present  a  factual  discussion  of  technical 
findings  and  accomplishments  for  the  quarter,  supported  by  technical 
descriptions,  mathematical  developments,  graphs,  and  illustrations  as  necessary. 
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4.  Site  Reports:  A  Site  Report  will  be  submitted  to  AFTAC  in  50  copies 
within  b0  days  after  the  completion  of  operations  at  each  site  visited. 

Each  report  will  be  identified  by  the  data  listed  in  paragraph  7  and  will 
include  the  notices  specified  in  paragraph  3.  Each  site  report  should  be 
comprehensive,  present ing  a  complete  and  factual  discussion  of  technical 
methods  and  findings  of  all  aspects  of  project  work  connected  with  the 
particular  site,  using  a  format  similar  to  that  of  the  site  and  analysis 
reports  published  under  Projects  VT/1139  and  VT/5051.  DD  Form  1473, 
Document  Control  Data  -  RAD  (reference  AFR  80-29)  will  be  included  in  each 
report.  The  appropriate  availability/ limitation  notice  for  use  on  the 
Form  1473  will  be  designated  by  AFTAC/TD-7. 

5.  Fir.'il  Report;  A  final  report,  Identified  by  the  data  listed  in 
paragraph  7  and  including  the  notices  specified  in  paragraph  8,  will  be 
submitted  to  AFTAC/TD-7  in  50  copies  within  60days  after  the  completion 
of  all  work.  The  report  will  present  a  concise  and  factual  discussion  of 
technical  findings  and  accomplishments  of  the  entire  project,  using  the 
site-report  Format.  DD  Form  1473,  Document  Control  Data  -  R&8  (reference 
AFR  80-29)  will  be  inoluded  in  each  report.  The  appropriate  availability/ 
limitation  notice  for  use  on  the  Form  1473  wi.,1  be  designated  by  AFTAC/TD-7 

6.  Special  Reports: 

a.  Special  reports  of  major  events  will  be  forwarded  by  telephone, 
telegraph,  or  separate  letter  as  they  occur  and  should  be  Included  in  the 
following  monthly  report.  Specific  items  are  to  include,  but  are  not 
restricted  to,  strikes  and  disasters,  program  delays,  technical  break¬ 
throughs,  major  decisions,  and  requirements  for  increase  in  funds. 

b.  Special  technical  reports  may  be  required  for  instrument  evalua¬ 
tions,  project  recommendations,  and  special  studies  when  it  is  more 
desirable  to  have  these  Items  reported  separately  from  the  periodic  and 
site  reports.  Specific  format,  content,  number  of  copies,  and  due  dates 
will  be  furnished  by  AFTAC  as  required. 

c.  All  seismograms  and  operating  logs,  including  pertinent 
information  concerning  time,  date,  type  of  instruments,  magnification, 
etc.,  will  be  provided  when  requested  by  the  AFTAC  project  officer. 

7.  Identification  Data;  All  monthly,  site,  and  final  reports  will  be 
identified  by  the  following  data: 

AFTAC  Project  No.  VELA  T/7703 

Project  Title:  Deep-Well  Array  Operations 

ARPA  Order  No.  624 

ARPA  Program  Code  No.  7F10 

Name  of  Contractor: 
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b.  All  site  and  final  reports  will  Include  the  following  notices  on 
the  first  or  title  page: 

This  research  was  supported  by. the 
Advanced  Research  Projects  Agency, 

Nuclear  Test  Detection  Office,  under 
Project  VELA-UNIFORM,  and  accomplished 
under  the  technical  direction  of  the 
Air  Force  Technical  Applications  Center 
under  Contract  No.  F33657-67-0-1224  » 

Qualified  users  may  request  copies  of  this  document  from: 

Defense  Documentation  Center 
Cameron  Station 
Alexandria,  Virginia  22314 
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